1/1 


L  HD-A13?  698 


UNCLASSIFIED 


IMPLEMENTATION  OF  NISITANI'S  BODV  FORCE  METHOD  FOR  THE 
SOLUTION  OF  NOTCH..  <U>  NATIONAL  AERONAUTICAL 
ESTABLISHMENT  OTTAWA  (ONTARIO)  W  E  FRAGA  ET  AL.  NOV  83 
NAE-AN-i?  NRC-22831  F/G  20/11 


NL 


ffi|\137690 


UNLIMITED 

UNCLASSIFIED 


Canada 


IMPLEMENTATION  OF  NISITANI’S 
BODY  FORCE  METHOD 
FOR  THE  SOLUTION  OF 
NOTCH  PROBLEMS 


by 

W.  E.  Fraga,  R.  L.  Hewitt 
National  Aeronautical  Establishment 


AERONAUTICAL  NOTE 


OTTAWA 
NOVEMBER  1983 


National  Research 
Council  Canada 


NAE-AN-17 
NRC  NO.  22831 


Conseil  national 
de  redherches  Canada 


Jw. 


NATIONAL  AERONAUTICAL  ESTABLISHMENT 
SCIENTIFIC  AND  TECHNICAL  PUBLICATIONS 


AERONAUTICAL  REPORTS: 

Aeronautical  Reports  (LR):  Scientific  and  technical  information  pertaining  to  aeronautics  considered 
important,  complete,  and  a  lasting  contribution  to  existing  knowledge. 

Mechanical  Engineering  Reports  (MS):  Scientific  and  technical  information  pertaining  to  investigations 
outside  aeronautics  considered  important,  complete,  and  a  lasting  contribution  to  existing  knowledge. 

AERONAUTICAL  NOTES  (AN):  Information  less  broad  in  scope  but  nevertheless  of  importance  as  a 
contribution  to  existing  knowledge. 

LABORATORY  TECHNICAL  REPORTS  (LTR):  Information  receiving  limited  distribution  because 
of  preliminary  data,  security  classification,  proprietary,  or  other  reasons. 

Details  on  the  availability  of  these  publications  may  be  obtained  from: 

Publications  Section, 

National  Research  Council  Canada, 

National  Aeronautical  Establishment, 

Bldg.  M-16,  Room  204, 

Montreal  Road, 

Ottawa,  Ontario 
K1A  0R6 

ETABLISSEMENT  AERONAUTIQUE  NATIONAL 
PUBLICATIONS  SCIENTIFIQUES  ET  TECHNIQUES 

RAPPORTS  D’AERONAUTIQUE 

Rapports  d’aeronautique  (LR):  Informations  scientifiques  et  techniques  touchant  l’aeronautique 
jugees  importantes,  completes  et  durables  en  termes  de  contribution  aux  connaissances  actuelles. 

Rapports  de  genie  mecanique  (MS).  Informations  scientifiques  et  techniques  sur  la  recherche  externe 
a  l’aeronautique  jugees  importantes,  completes  et  durables  en  termes  de  contribution  aux  connais¬ 
sances  actuelles. 

CAHIERS  D’AERONAUTIQUE  (AN):  Informations  de  moindre  portee  mais  importantes  en  termes 
d’accroissement  des  connaissances. 

RAPPORTS  TECHNIQUES  DE  LABOR ATOIRE  (LTR):  Informations  peu  disseminees  pour  des 
raisons  d’usage  secret,  de  droit  de  propriete  ou  autres  ou  parce  qu’elles  constituent  des  donnees 
preliminaires. 

Les  publications  ci-dessus  peuvent  etre  obtenues  a  Fadresse  suivante: 

Section  des  publications 

Conseil  national  de  recherches  Canada 

Etablissement  aeronautique  national 

Im.  M-16,  piece  204 

Chemin  de  Montreal 

Ottawa  (Ontario) 

K1A  0R6 


UNLIMITED 

UNCLASSIFIED 


IMPLEMENTATION  OF  NISITANI’S  BODY  FORCE  METHOD 
FOR  THE  SOLUTION  OF  NOTCH  PROBLEMS 


APPLICATION  DE  LA  METHODE  DES  FORCES  DE  NISITANI  POUR 
RESOUDRE  DES  PROBLEMES  D’ENTAILLES  ET  DE  FISSURES 


by /par 

W.E.  Fraga*,  R.L.  Hewitt 
National  Aeronautical  Establishment 

*  Department  of  Mathematics,  University  of  Alberta 


AERONAUTICAL  NOTE 

OTTAWA  NAE- AN-17 

NOVEMBER  1983  NRC  NO.  22831 


W.  Wallace,  Head/Chef 

Structures  and  Materials  Laboratory/ 

Labors toire  des  structures  et  materiaux 


G.M.  Lindberg 
Director/Directeur 


SUMMARY 


'  -y  An  outline  of  Nisitani’s  body  force  method  is  presented  and 
details  of  the  numerical  implementation  provided.  Full  details  are  provided 
for  the  solution  of  several  problems  including  single  cracks  or  single  elliptical 
notches  in  semi-infinite  plates  as  well  as  the  problem  of  an  embedded  crack 
near  a  semi-elliptical  notch  in  a  semi-infinite  plate. 


RESUME 


Une  description  de  la  methode  des  forces  de  Nisitani  avec  des 
exemples  numeriques  est  presentee  dans  cette  note.  Les  exemples  de 
resolution  de  problemes  sont  tres  detailles  et  comprennent  notamment  des 
problemes  de  fissures  uniques  ou  d’entailles  elliptiques  uniques  dans  des 
plaques  semi-infinies  ainsi  que  le  probleme  d’une  fissure  situee  a  proximite 
d’une  entaille  semi-elliptique  dans  une  plaque  semi-infinie. 
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IMPLEMENTATION  OF  NISITANrS  BODY  FORCE  METHOD 
FOR  THE  SOLUTION  OF  NOTCH  PROBLEMS 


1.0  INTRODUCTION 

Efficient  use  of  fracture  mechanics  methods  for  damage  tolerance  assessments  depends  upon 
knowledge  of  the  stress  intensity  factor  for  the  flaw  considered.  While  these  are  available  for  some 
simple  crack  geometries,  it  is  often  necessary  either  to  resort  to  some  approximation  or  employ  a 
complex  finite  element  program  to  obtain  the  information.  Alternatively,  a  boundary  integral 
approach  may  be  used  which  may  in  some  cases  be  more  efficient. 

One  such  approach  is  the  body  force  method  proposed  by  Nisitani!  .  This  is  a  relatively 
straightforward  technique  in  principle  which  makes  use  of  the  stress  field  derived  from  point  forces 
acting  in  an  infinite  plate  or  body.  Boundary  conditions  are  satisfied  by  applying  body  forces 
(continuously  embedded  point  forces)  along  the  edges  of  notches  or  cracks  (the  density  of  these  body 
forces  being  determined  from  the  boundary  conditions).  The  method  is  applicable  to  both  two  and 
three  dimensional  problems. 

In  attempting  to  use  this  method  following  the  work  of  Nisitani,  it  was  found  that  certain 
details  of  the  technique  were  absent  or  available  only  in  the  Japanese  literature,  and  it  therefore  seems 
appropriate  that,  having  implemented  the  method  at  NRC,  the  details  of  the  implementation  be 
reported.  Thus  for  the  particular  problems  studied,  all  the  relevant  information  from  various  papers  of 
Nisitani  have  been  brought  together  and  all  the  details  necessary  for  implementation  are  given  in 
this  report.  An  additional  report!2*  is  available  on  request  which  contains  listings  of  the  computer 
programs  developed  with  appropriate  documentation. 


2.0  PRINCIPLES  OF  BODY  FORCE  METHOD 

In  order  to  make  this  report  more  complete,  a  brief  review  of  the  principles  of  the  body 
force  method  will  be  presented.  The  review  follows  that  of  Nisitani!1!  but  various  points  which  were 
found  difficult  to  follow  have  been  amplified  while  others  which  were  considered  to  be  unnecessary 
for  a  general  understanding  of  the  method  have  been  omitted. 

2.1  Tension  of  Infinite  Plate  with  an  Elliptic  Hole 

Consider  the  problem  of  an  elliptical  hole  in  an  infinite  plate  as  shown  in  Figure  1.  If  a 
stress,  o“,  is  applied  at  infinity,  the  x  and  y  components  of  the  displacement,  u  and  v,  at  an 
arbitrary  point  ({,  r?)  on  the  edge  of  the  hole  are  independent  of  the  Poisson’s  ratio,  v,  and  are  given 
byl3l : 


u  -  —  (1  +  2b/a)  * 


where  E  is  the  Young’s  modulus  and  b  and  a  are  the  semi-major  and  semi-minor  widths  respectively 
of  the  ellipse. 


V-V- ' N'.  V.  ■ 


Now  consider  an  elliptical  plate  of  the  same  size  and  shape  as  the  hole.  For  the 
displacements  on  the  edge  of  this  plate  to  be  the  same  as  for  the  edge  of  the  hole,  the  strain 
components  at  an  arbitrary  point  on  the  plate  must  be 


du  _  °x 

ai  ”  ~e 


ex  =  TZ  =  -(l  +  2b/a) 


7vy  =  0 


where  e  and  y  are  the  normal  and  shear  strains  respectively.  From  Hooke’s  law,  the  normal  and  shear 
stress,  a  and  r,  corresponding  to  these  strains  are 


- -  (ex  +  vey)  =  - -  (1  +  2b/a  -  v) 

1-v2  1  -v2 


- -  (ey  +  vex) - (1  -  *  (1  +  2b/a)) 

1  -  v2  1-v2 


T\y  ~  0 


The  elliptical  plate  having  these  surface  stresses  can  be  inserted  exactly  into  the  infinite  plate  with 
the  elliptical  hole  subjected  to  a  stress  o“,  producing  an  infinite  plate  with  no  hole.  Therefore  the 
problem  of  an  infinite  plate  with  a  hole  is  equivalent  to  the  problem  of  an  infinite  plate  with  no  hole, 
but  having  the  stresses  given  by  Equations  (3)  applied  along  the  imaginary  boundary  of  the  ellipse  as 
shown  in  Figure  2.  These  stresses  are  obtained  from  a  series  of  continuously  embedded  point  forces 
along  the  boundary. 

The  densities,  px  and  py,  of  these  body  forces  in  the  x  and  y  directions  can  be  obtained 
from  Equations  (3)  as 


px  =  - -  (1  +  2b/a  -  v) 

1-v2 


py - ;  (1  -  *  (1  +  2b/a)) 

1-v2 


where  ‘unit  length’  is  measured  in  y  direction  for  px  and  in  x  direction  for  pv 


Since  v  does  not  affect  the  end  results  according  to  Nisitanil 1  * ,  it  may  be  set  to  zero  and  the  body 
forces  written 


px  =  o“  (1  +  2b/a) 


(5a) 


and 


(5b) 


Thus  for  an  elliptical  hole  in  an  infinite  plate  the  densities  of  body  force  are  constant  along  the 
boundary  of  the  hole.  In  general,  however,  they  will  not  be  constant  and  a  numerical  technique  must 
be  used  to  determine  them. 


2.2  Tension  of  Infinite  Plate  with  a  Crack 

The  case  of  an  infinite  plate  with  a  crack  under  tension  is  shown  in  Figure  3.  A  crack  can 
be  considered  as  an  extremely  slender  elliptical  hole.  As  b/a  tends  to  infinity,  the  hole  becomes 
a  crack  and  the  points  at  which  body  forces  are  applied  come  closer  to  each  other.  Thus  the 
fundamental  stress  field  is  obtained  from  a  pair  of  point  forces  equal  in  magnitude  and  opposite  in 
direction. 


To  find  the  stress  at  an  arbitrary  point  (x,  y)  due  to  a  pair  of  point  forces  separated  by  a 
distance  5  from  each  other,  consider  Figure  4.  The  stress  at  (x,y)  in  the  x  direction  due  to  the  body 
forces  X  acting  on  the  element  d?j  may  be  written 


o 


X 


x  a?  ($,T?,X, y) 


X=1 


(6) 


where  a*  is  the  stress  in  the  x  direction  at  (x,y)  due  to  a  unit  body  force  in  the  x-direction  at 
({,17).  This  may  be  written 


a 


X 


•  5 


3o? 

a* 


(X«) 


(7) 


But  the  point  force  X  can  be  written 


and  for  the  very  slender  ellipse  (i.e.  small  {) 


6  -  2* 


Therefore,  substituting  Equations  (8)  and  (9)  into  (7) 


But  for  an  ellipse 


and  from  Equation  (5a) 


°x  *  (P*dT?  2f) 


dox  =  lim  J(pxdT?) 


3oxx 

~rr  Px  dr? 

x  =  i  t — o 

t=o 


£  +  £.i 

a2  b2 


Px  =  ffx  ^  +  2b/a> 


Thus  substituting  Equations  (5a)  and  (11)  in  (10) 


dax  =  •  4a"  V  b2  -  tj2  d rj 


From  the  fundamental  stress  field  for  a  point  force  acting  in  an  infinite  plate  (see  Appendix  1) 


-(x-t)(3(x-£)2  +  (y-r?)2) 


4  ir  ((x-f)2  +  (y-r?)2) 


therefore 


(14 


ooc 


-  3x4  +  6x2  (y  - 17  )2  +  (y  - 1?)4 


x  =  i 

E  =  0 


4tt(x2  +  (y-1?)2  )3 


To  determine  the  stress  at  an  arbitrary  point  (0,y)  on  the  y  axis  (i.e.  on  the  crack  line) 
Equation  (12)  must  be  integrated  and  added  to  the  uniform  stress  at  infinity  0“ ,  i.e. 


crx(y)  =  lira 

x-»0 


f  (4o:v^v)d, 

'  4tt(x2  +  (y-17)2  )3  ' 

_r  b 


+  a: 


(i5: 


When  1  y  |  >  b,  i.e.  beyond  the  crack,  x  may  be  set  to  zero  directly  and  Equation  (15)  may 

be  written 


Vb2-r)2 

(y-17)2 


drj  +  a“ 


iyl 


sin-1 


\/y2  -  b2 


yr?~b2 
b(y  - 17) 


Iyiq: 


0“  Vffb 

— -  —  for  r  =  |  y  |  -  b  <  b 

v/2^ 


(16 


or  ox(y)  =  _  (17 

V  2ttt 


where  Kj  is  the  mode  1  stress  intensity  factor 

For  |y  |  <  b  and  x  =  0,  the  integral  contains  a  singularity  at  y  =  tj,  so  principal  values 
must  be  taken  to  evaluate  it.  Denoting  the  integrand  as  /(ry).  Equation  (15)  may  be  written 


a 


& 


I 


ox(y)  =  lim  I  /(r?)  drj  +  a 

x-0  J  * 

-b 


x-0  J 

_-b 


v 

+  lim  |  /(r?)  dr?  +  o“ 

x  -*0  J 


Hv)  di?  +  /(T?)dr? 


The  first  two  integrals  contain  no  singularities  and  may  be  integrated  directly  by  setting  x  =  0.  The 
result isl1!  for  |e|  <  b: 

-  y-e  b  “j 

Urn  J  /(Tj)dTj  +  j  /(i?)  dr?  I  =  [F(y-e)  -  F(-b)  +  F(b)  -  F(y+e)] 

L-b  y+e  J 


JL  1 

it 


2  Vb2  -  y2 


in  which 


_°o 

■if 


>/b2  -  t?2 

(y  -  t?)2 


^  f^Eil .  sin-i  (i) 

■n  I  y-v  \b / 


y  .  (b2  -  yr?)  +  V(b2  -  y2)  (b2  - 

+ - In  - 

v/JTTJT  b(y-,) 


~j] 


The  third  integral  contains  a  singularity  when  x  -  0,  so  it  is  not  permissible  to  set  x  =  0 
before  performing  the  integration.  The  integration  must  be  performed  first  before  the  limit  is  taken. 

For  e  <  b,  Vb2-tj2  may  be  approximated  by  Vb2  -  y2  which  may  then  be  taken  outside  the 
integral.  Then  setting  t  -  y  -  rj 


-7  - 


I 


•3 


lim  /  /(t?)  dt?  —  lira 

x-»0  J  x-»0 


2  _  „2 


-  y 


+  6x2t2  +  t4  ~l 
- dt 

c2+t2)3  J 


=  lim  I  — 

x"°  L 


'pTy1  -i-  +  -Jg5- 

(x2+e2  (x2+e2)2 


—  V  b2  -  y2 
ire 


Substituting  Equations  (19)  and  (21)  into  (18)  yields  the  result 


ax(y)  =  0  for 


I  y  I  <  b 


Equations  (22)  and  (17)  show  that  the  problem  of  an  infinite  plate  with  a  crack  can  be 
reduced  to  that  of  an  infinite  plate  without  a  crack  but  subjected  to  body  forces  along  the  imaginary 
crack  line.  Hie  stresses  at  an  arbitrary  point  are  obtained  by  integrating  the  stress  field  for  a  pair  of 
point  forces. 

It  is  also  clear  from  examination  of  Equation  (20)  that  the  same  result  would  be  obtained 
by  ignoring  the  singularity,  and  in  fact  no  special  consideration  is  necessary  when  integrating  these 
expressions  analytically. 

Hie  strength,  dR,  of  a  pair  of  point  forces  has  been  expressed  as  the  product  of  the 
magnitude  of  the  forces,  X,  and  the  distance  between  them,  6, 


dR  *  X5 


or,  from  Equations  (8)  and  (9) 


dR  =  pxdr?2$  I 


Hien  the  density  of  a  pair  of  point  forces,  dR/drj,  is  obtained  as 


dR 

ST  ’  <p«2t) 


which  for  a  crack  in  an  infinite  plate  may  be  written  using  Equations  (5a)  and  (11)  as 


-j-  ■  4 o“  Vb 2  -  t?2 

dr?  x 


Ni-i 

..  a 


However,  for  convenience,  Nisitani^ 1 ' ,  has  defined  a  quantity,  7,  which  for  the  crack  in 
the  infinite  plate  is  given  by 


1 

4\/b2  - 172 


dR 

drj 


(24) 


and  which  he  calls  the  density  of  a  pair  of  body  forces. 

(The  reason  for  using  this  quantity,  rather  than  dR/d rj  directly,  is  that  the  variation  in  7  along  the 
crack  line  for  general  problems  will  then  be  small,  and  so  setting  7  =  constant  within  an  interval  will 
result  in  an  accurate  solution). 

Comparing  Equations  (23)  and  (24),  for  the  crack  in  an  infinite  plate,  7  is  constant  and  equal  to  o“. 
Thus  from  Equation  (17) 


Kj  =  y\/lrb 


(25) 


For  other  crack  problems,  the  densities  of  pairs  of  body  forces  7  will  not  be  constant,  but  will  vary 
along  the  crack.  For  these  problems  a  numerical  method  must  be  used  to  determine  the  densities  of 
body  forces  at  the  tips  of  the  cracks  yA  and  7B.  From  these  the  stress  intensity  factors  K1A  and  K,B 
at  each  tip  can  be  calculated  according  to  Nisitani* 1 '  as 


k\a  =  7A\/ffb 


(26) 


^1B  =  7b  V^b 


3.0  NUMERICAL  TECHNIQUE 
3.1  Crack  in  a  Semi-infinite  Plate 
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3.1.1  Formulation 


Consider  the  case  of  a  crack  in  a  semi-infinite  plate  under  tension  as  shown  in  Figure  5. 
Whereas  in  an  infinite  plate  the  densities  of  pairs  of  body  forces,  7,  were  constant  along  the  crack, 
here  the  densities  vary  continuously. 


If  the  crack  is  divided  into  MM  equal  intervals  and  7  is  set  constant  within  each  interval,  a 
stepped  distribution  will  result  which  will  be  accurate  for  large  values  of  MM.  The  derivative  of  the 
stress  field  due  to  a  point  force  in  a  semi-infinite  plate  is  used  to  determine  the  stresses  at  the  mid¬ 
point  of  each  interval  due  to  the  body  forces  and  the  boundary  conditions  that  the  mid-point  of  each 
interval  be  free  of  stress  are  applied.  The  densities  of  the  pair  of  body  forces  in  each  interval  can  then 
be  determined. 


The  stress  intensity  factors  at  both  tips  of  the  crack  can  be  calculated  using  only  the 
densities  7t  and  7mm  at  the  ends  of  the  crack.  The  stresses  at  an  arbitrary  point  on  the  plate  can  be 
found  by  taking  a  linear  combination  of  the  body  force  densities. 
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3.1.2  Solution  Sequence 

(1)  The  crack  AB  is  divided  into  MM  intervals,  with  the  end  points  of  the  N-th  intervals  having 
values  of  i?  of 


qN_j  =  e  -  b  +  2(N  -  l)b/MM 


%  =  e  -  b  +  2Nb/MM 


where  b  is  the  half  length  of  the  crack  and  e  the  distance  from  the  centre  of  the  crack 
to  the  edge  of  the  plate  as  shown  in  Figure  5.  The  mid-point  of  the  M-th  interval  has  a 
co-ordinate,  y  of 


yM  =  e-b  +  (2M  -  l)b/MM 


(2)  The  influence  coefficient  o^,  which  is  the  stress  at  the  mid -point  of  the  M-th  interval  due 
to  a  pair  of  body  forces  on  the  N-th  interval  having  a  density  7N  =  1 ,  is  given  by 


77n  = 

xM 


(  Ha* 

J 


N-th 
interval 


4Vb2  -  (tj  -  e)2  dr? 


XM 
x  =  t  =0 


Using  the  fundamental  stress  field  for  a  point  force  acting  in  a  semi-infinite  plate 
[Appendix  1]  this  may  be  written 


1  12y  12y2 

(y  +  q)2  (y  +  r?)3  (y  +  r?)4 


—  •  \/b2  -  (rj  -  e)2  drj 

it 


This  integral  can  be  solved  either  analytically  (see  3.1.3  and  3.1.5)  or  numerically,  with  the 
singularity  treated  in  the  same  manner  as  for  the  infinite  plate. 

(3)  The  boundary  conditions  are  applied  to  make  the  mid -point  of  each  interval  free  of 
stress,  i.e. 


MM 

N=  1 


for  M  =  1,2  . .  .MM 


This  results  in  a  system  of  MM  linear  equations  in  MM  unknowns,  ?N ,  which  can  be  solved 
using  a  matrix  inversion  routine. 


>  j  a  ,s  %  v 


(4)  The  stress  intensity  factors  K1A  and  K1B  are  then  calculated  from  7,  and  7MM  as 

kia  =  and  K1B  =  7mm 

(5)  The  stress  at  an  arbitrary  point  P  on  the  plate  can  be  calculated  from  the  influence 
coefficients  for  that  point  and  the  body  force  densities,  i.e. 


axP  = 


£ 


Tn  °x?  +  °x 


(6)  For  finite  values  of  MM  there  will  be  errors  in  K1A  and  K1B  because  of  the  approximation 
of  a  stepped  distribution.  The  true  values  of  K1A  and  K1B  correspond  to  MM  =  Since 
Nisitani  has  shown  the  error  to  be  proportional  to  1/MM  these  values  can  be  obtained  by 
extrapolation  from  two  values  of  MM  such  as  12  and  24. 

3.1.3  Analytic  Integration 

The  integral  in  3.1.2  in  the  calculation  of  was  performed  analytically.  The  result  is 
shown  in  Appendix  2.  Although  the  integral  involves  a  singular  term  in  the  case  where  M  =  N,  no 
special  consideration  is  required  as  illustrated  earlier. 

3.1.4  Results  (Analytic  Integration) 

Table  I  compares  the  correction  factors,  F1A,  obtained  in  the  present  investigation  with 
those  of  Nisitani^ 1 '  for  various  values  of  the  number  of  intervals  MM,  for  a  crack  with  a  b/e  ratio 
of  0.5. 


TABLE  I 


CRACK  IN  SEMI-INFINITE  PLATE 


12  1.090367  1.09033  b/e  =  0.5 

16  1.090622  1.09063 

24  1.090919  1.09091  K 

32  1.091024  1.09104  F1A  =  - — 

48  1.091243  1.09115  o“\/»rb 

The  agreement  is  very  good  and  the  minor  discrepancies  can  be  attributed  to  computer  differences. 

The  correction  factor  F1A  is  plotted  against  1/MM  in  Figure  6.  The  error  in  Fj  A  is  almost 
linear  with  respect  to  1/MM,  with  the  value  of  F,  A  tending  to  1.09149  as  MM  tends  to  infinity.  This 
compares  with  the  value  of  1 .0914  obtained  by  Nisitani1 1 1 . 


v  v 


v  A  --.  '. 
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3.1.5  Numerical  Integration 

For  most  problems,  the  integral  for  the  influence  coefficients  cannot  be  integrated 
analytically  as  in  3.1.3  and  the  integration  must  be  performed  numerically.  To  check  on  suitable 
parameters  for  these  integrations,  the  integral  of  3.1.2  used  to  evaluate  °xM  was  also  solved  using  a 
numerical  integration  routine. 

When  evaluating  the  influence  coefficients  for  M  =  N,  i.e.  the  stresses  at  interval  M  due  to 
unit  body  forces  in  the  same  interval,  there  is  a  singularity  which  must  be  taken  into  consideration 
when  using  numerical  integration.  The  integrals  either  side  of  the  singularity  are  evaluated  in  the 
normal  manner  and  the  integrand  is  approximated  in  the  area  of  the  singularity  as  follows: 


y+e  r 


F(tj) 


-7  / 


y-e  «- 


_1 _  12y  12y2 


(y  ~n)2  (y+*?)2  (y+n)3  (y  +  n)4 


\/b2  -  (17-  e)2  dtj 


For  small  e  this  becomes 


y  +  e 


1  r  i  _ 

F(t7)  ~  —  /  - Vb2  -  (i?-  e)2  drj 

*  J  (y-i?)2 


y-e 


Also  for  small  e,  tj  ~  y 


y  +  e 


F(t?)  -  -  Vb2-  (y-e)2  f  - - —  dr) 

*  J  (y-r?)2 

y-e 


2_ 

ire 


Vb2  -  (y  -  e)2 


Thus  the  total  integral  over  an  interval  [r?N  _  j ,  tjn  ]  which  contains  a  singularity  can  be  evaluated 
numerically  as  the  sum  of  three  integrals  over  three  sub-intervals,  the  outer  two  being  performed 
using  a  standard  numerical  integration  routine,  while  the  middle  one  is  a  simple  analytic  expression. 

3.1.6  Influence  of  Integration  Parameters  and  Computer  Precision 

The  influence  coefficients  which  involved  the  singularities,  ,  were  found  to  be  affected 
by  the  choice  of  e  in  the  numerical  integration.  The  error  was  linear  with  respect  to  e,  however,  as 
shown  in  Table  II,  so  extrapolation  of  the  true  value  was  straightforward. 

The  integration  parameter  e  was  always  calculated  as  a  linear  function  of  the  fraction  of 
the  interval  (FRAC)  which  is  to  be  evaluated  analytically  around  the  singularity.  In  this  way,  as  the 
intervals  become  smaller  (as  MM  is  increased)  the  value  of  e  shrinks  accordingly  without  any  change 
to  FRAC. 


TABLE  II 


INFLUENCE  OF  e  ON  ERROR  IN  a^1 


FRAC 

e 

ERROR  IN  al\ 

0.020 

0.0020 

0.00748 

0.018 

0.0018 

0.00673 

0.016 

0.0016 

0.00599 

0.014 

0.0014 

0.00525 

0.012 

0.0012 

0.00451 

0.010 

0.0010 

0.00377 

0.008 

0.0008 

0.00304 

0.006 

0.0006 

0.00230 

0.004 

0.0004 

0.00160 

0.002 

0.0002 

0.00095 

Results  for  crack  in  infinite  plate  in  the  interval  [1.0, 1.2] 
e  =  2.0;  b  =  1.0;  y  =  1.1;  MM  =  10;  e  =  FRAC  •  b/MM 


The  correction  factors  F1A  and  F1B  where 


_Kl_ 

o“V^b 


and  stress  intensity  factors  were  also  found  to  vary  linearly  with  respect  to  e,  so  for  a  given  MM, 
three  or  four  values  of  e  were  used  in  the  numerical  integration  and  the  true  values  were  extrapolated 
using  linear  regression. 

Hie  number  of  digits  of  precision  specified  for  real  program  variables  in  the  computer 
calculations  was  found  to  have  a  great  effect  on  the  results.  For  the  NRC  IBM  360  system,  a  ‘single 
precision’  real  variable  has  six  digits  of  precision,  while  a  ‘double  precision’  variable  has  fourteen.  For 
these  calculations  it  was  necessary  to  specify  double  precision  for  all  real  variables. 

Hie  primary  reason  for  this  is  that  for  small  values  of  e,  several  digits  of  precision  are  lost 
in  the  integration  over  intervals  containing  the  singularities.  When  evaluating  the  integral  over  the 
three  sub-intervals,  the  outer  ones  produced  large  positive  numbers,  while  the  centre  produced  a 
large  negative  number.  When  the  three  were  added,  the  result  was  a  number  close  to  zero.  Thus  in  the 
addition  three  or  more  digits  of  precision  were  lost. 


3.1.7  Influence  of  Number  of  Intervals 

Figure  7  shows  that  the  error  in  the  correction  factor  or  stress  intensity  factor  was  linear 
with  respect  to  1/MM,  and  true  values  were  therefore  obtained  by  extrapolation  to  MM  -  °°. 


3.1.8  Results  (Numerical  Integration) 


Table  III  compares  the  results  obtained  by  numerical  integration  with  those  obtained 
analytically  for  various  values  of  MM.  The  agreement  is  good  and  discrepancies  can  be  attributed  to 
approximation  errors  in  the  integration  routine. 


TABLE  III 

COMPARISON  OF  F,A  OBTAINED  BY  NUMERICAL  AND  ANALYTIC 
INTEGRATIONS  FOR  CRACK  IN  SEMI-INFINITE  PLATE 


MM 

FRAC 

(EXTRAPOLATED) 

0.00 

ANALYTIC 

INTEGRATION 

0.03 

0.02 

0.01 

12 

1.096155 

1.094206 

1.092265 

1.090319 

1.090367 

16 

1.096337 

1.094411 

1.092468 

1.090536 

1.090622 

24 

1.096515 

1.094635 

1.092756 

1.090876 

1.090919 

32 

1.096554 

1.094688 

1.092795 

1.090920 

1.091024 

48 

1.096608 

1.094777 

1.092938 

1.091104 

1.091243 

b/e  =  0.5  e  =  FRAC  •  b/MM 


The  correction  factor  F1A  as  e  tends  to  zero  is  plotted  against  1/MM  in  Figure  7.  Once 
again,  the  error  in  F1A  is  nearly  linear  against  1/MM  and  the  value  tends  to  1.091352  compared 
with  a  value  of  1 .091490  using  analytic  integration  and  1 .0914  obtained  by  Nisitani' 1  • . 

3.2  Edge  Crack  in  a  Semi-infinite  Plate 
3.2.1  Formulation 

The  case  of  a  crack  at  the  edge  of  a  plate  as  in  Figure  8  is  similar  to  the  case  of  an  embedded 
crack.  The  edge  of  the  plate  is  placed  at  the  centre  of  the  crack,  so  that  the  distance,  e,  from  the  edge 
to  the  centre  of  the  crack  is  zero. 


3.2.2  Results  (Numerical  Integration) 

Table  IV  shows  the  correction  factor  F1B  for  various  values  of  the  number  of  intervals, 
MM,  and  the  integration  parameter,  FRAC,  together  with  the  extrapolated  values  for  FRAC  =  0. 

The  correction  factor  F1B  is  plotted  against  1/MM  in  Figure  9.  The  error  in  F1B  is  almost 
linear  in  1/MM  with  the  value  tending  to  1.1213  as  MM  tends  to  infinity.  This  compares  with  a  value 
of  1.121  given  by  Nisitanil1' . 


VSN 


TABLE  IV 

COMPARISON  OF  F,b  FOR  EDGE  CRACK  IN  SEMI-INFINITE  PLATE 
AT  VARIOUS  VALUES  OF  MM 


FRAC 

8 

12 

16 

24 

0.03 

1.096820 

1.098911 

1.099955 

1.101096 

0.02 

1.101942 

1.104243 

1.105375 

1.106618 

0.01 

1.107365 

1.109990 

1.111269 

1.112726 

0.00 

(EXTRAPOLATED) 

1.112587 

1.115460 

1.116847 

1.118443 

e  =  FRAC  •  b/2MM 


Fib  = 


a”  y/irb 


3.3  Elliptical  Hole  in  an  Infinite  Plate 


The  numerical  solution  of  problems  involving  holes  is  a  little  more  complex  than  for  cracks. 
It  was  therefore  decided  to  first  try  the  simple  problem  of  an  elliptical  hole  in  an  infinite  plate  since 
the  body  force  densities  are  known  and  constant  in  this  case,  and  this  allows  an  intermediate  check  of 
the  program. 

3.3.1  Formulation 

The  elliptical  hole  shown  in  Figure  1  is  divided  into  MM  equal  intervals  and  the  body  force 
densities,  px  and  py  are  set  constant  within  each  interval.  The  influence  coefficients  are  calculated 
from  the  stress  field  for  a  point  force  in  an  infinite  plate  and  the  boundary  conditions  that  the  mid¬ 
point  of  each  interval  be  free  from  stress  are  applied  to  determine  the  body  force  densities.  The  stress 
concentration  factors  and  the  stresses  at  an  arbitrary  point  on  the  plate  can  be  found  from  a  linear 
combination  of  the  body  force  densities. 

3.3.2  Solution  Sequence 

(1)  The  ellipse  is  divided  into  MM  equal  intervals,  using  the  elliptic  parameter  p  used  in  the 
equations  of  an  ellipse: 

£  =  a  cos  p  ;  r?  =  b  sin  p 
The  values  of  <f>  at  the  end  points  of  the  N-th  interval  are 

PN-i  =  -ir/2  +  2(N-1)*/MM 

and  <pN  *  -ir/2  +  2Nir/MM 


Hie  value  of  <p  at  the  mid-point  of  the  N-th  interval  is 


*  -  -  ff/2  +  (2N- 

The  x  and  y  values  of  the  mid-point  of  the  N-th  interval  are 


x  =  acos(-7r/2  +  (2N-1)*/MM) 
y  =  b  sin  (-  jr/2  +  (2N-1)tt/MM) 

(2)  The  influence  coefficients,  ;  o™ ;  o*  *;  o™ ;  r*yNM ;  ,  which  are  the  stresses  at 

the  mid-point  of  the  M-th  interval  due  to  a  body  force  on  the  N-th  interval  having  densities 
px  *  1  or  py  =  1  are  given  in  Appendix  3. 

In  the  program  listed  in  Reference  2,  the  influence  coefficients  for  each  point  on  the  ellipse 
have  been  calculated.  However,  because  of  symmetry  about  the  x  and  y  axes  the  influence 
coefficients  need  only  be  calculated  for  one  quarter  of  the  ellipse,  although  for  that  quarter, 
all  the  influence  coefficients  must  be  calculated,  i.e.  the  effects  ‘of  body*  forces  around 
the  total  ellipse  must  be  calculated. 

(3)  Hie  boundary  conditions  to  make  the  mid-point  of  each  interval  free  from  stresses  are 
applied;  i.e. 


MM 

Px N  (°xM  COS  0  +  T?yM  •“»  0) 

N  =  1 


MM 


-  0 


and 


MM 

Px  N(°yXMN  8in  6  +  rxyM  COS  0) 

N  » i 


MM 


.  /_YN 


where  pxN,  pyN  are  the  body  force  densities  of  the  N-th  interval,  o“  is  the  stress  at 
infinity  and  0M  is  the  angle  between  the  x  axis  and  the  normal  to  the  ellipse  at  the  mid¬ 
point  of  the  M-th  interval;  i.e. 


where  ^  is  the  value  of  <p  at  the  mid-point  of  the  M-th  interval. 

This  results  in  a  system  of  2MM  linear  equations  in  2MM  unknowns  pxN,  PyN- 

(4)  The  stresses  at  any  point  P  in  the  plate  are  calculated  from  the  influence  coefficients  at 
that  point  and  the  body  force  densities;  i.e. 

MM 

••  -  Y  (<4K  *  C  p,h)  -  •; 

N=  1 


MM 

„  _  \  '  /  XN  „  ,  „YN  „  \ 
°y  2^  (CTyP  PxN  +  °yP  PyN/ 

N=  1 


MM 

’  Y  '«  +  T™’  '*) 

N=1 


(5)  The  stress  concentration  factors  KJA,  Kj.b  are  then  calculated: 

KrA  =  ox  (at  A)/o” 

“tb  "  B)/a“ 

(6)  Because  the  body  force  densities  are  constant  around  the  ellipse,  there  will  be  no  error  due 
to  the  choice  of  the  number  of  intervals  MM.  Therefore  the  value  of  the  stress  concentration 
factor  will  be  the  same  regardless  of  the  number  of  intervals. 

3.3.3  Numerical  Integration 

In  the  numerical  integration,  the  stress  field  equations  have  a  singularity  at  y  =  rj  and 
x  ■  £,  i.e.  when  the  stresses  are  evaluated  at  the  mid-point  (x,y)  of  an  interval  due  to  a  point  force 
X  or  Y  on  the  same  interval.  Integration  is  straightforward  over  all  other  intervals. 

In  order  to  integrate  over  the  singularity,  principal  values  of  the  integral  are  taken  by 
splitting  the  integral  in  three  disjoint  subintervals 

l*N-I'*MID“eJ’’  ^MID"e'^MID  +  e1*  ^MID  +  e*^Nl 


where 


x  =  a  cos  <pMID;  y  =  bsin^M,D 
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The  integrals  F !  (v?)  and  F3  (i p)  over  the  1st  and  3rd  subintervals  are  straightforward  because 
the  singularity  occurs  in  F2  over  the  2nd  subinterval. 

If  e  is  small  relative  to  the  length  of  the  interval  [ipN  _  j ,  >PN  ] ,  F2  (<p)  can  be  approximated 
by  a  function  which  contains  no  singularity  (see  Appendix  4). 

For  this  same  M-th  interval  when  x  =  £  and  y  =  17,  the  stresses  A ax ,  AaJ,  Aa* ,  AaJ, 
At*,,  At*,  which  are  the  stresses  at  a  point  infinitesimally  close  to  the  ellipse  subjected  to  a  body 
force  must  be  added  to  o*^,  °yM>  ayM»  rxyM»  TxyM  respectively. 

It  can  be  shown  (Appendix  5)  that  these  additional  stresses  are 


- (5  +  4  cos  20  -  cos  40) 

16 


Aa*  =  — 

y  „  16 
"x-1 


(1  -  cos  40) 


- (2  sin  20  -  sin  40) 

16  ’ 


Aax  I  =— (1-COS40) 


Aa*  = - (5-4  cos  20  -  cos  40) 

y  16 
Py-1 


- (2  sin  20  +  sin  40 ) 

16 


where  0  is  the  angle  between  the  x  axis  and  the  normal  to  the  ellipse  at  the  mid-point  of  the  M-th 
interval. 

3.3.4  Results 

Table  V  shows  the  body  force  densities  px  and  py  for  each  interval  of  a  circular  ellipse  in 
an  infinite  plate. 

The  results  are  identical  to  those  predicted  by  the  theory  in  Section  2.1.  The  body  force 
densities  are  constant  around  the  ellipse  and  agree  with  equations 


px  =  0“  (1  +  2  b/a)  =  3o~  =  3 


Py  =  ~  °x  ~  ^ 


TABLE  V 


BODY  FORCE  DENSITIES  FOR  CIRCLE  IN  INFINITE  PLATE 


Interval 

Px 

Py 

1 

3.000000 

-1.000000 

2 

3.000000 

-l.oooono 

3 

3.000000 

-1.000000 

4 

3.000000 

-1.000000 

5 

3.000000 

-1.000000 

6 

3.000000 

-1.000000 

7 

3.000000 

-1.000000 

8 

3.000000 

-1.000000 

MM  =  8;  FRAC  =  0.02;  e  =  FRAC  •  tt/MM 
a  =  1.0;  b  =  1.0;  a”  =  1.0 


3.4  Elliptic  and  Semi-elliptic  Hole  in  a  Semi-infinite  Plate 

3.4.1  Formulation 

For  the  case  of  an  ellipse  or  a  semi-ellipse  in  a  semi-infinite  plate.  Figure  10,  the  densities 
of  body  force  px  and  py  are  no  longer  constant  as  in  the  case  of  the  ellipse  in  an  infinite  plate. 
Instead  the  densities  are  continuously  varying  around  the  edge  of  the  hole. 

If  the  elliptic  hole  is  partitioned  into  MM  intervals  and  px,  py  are  set  constant  within  each 
interval,  the  continuously  varying  densities  will  be  replaced  by  a  stepped  distribution  which  will  be 
accurate  for  large  values  of  MM. 

Tlie  stress  field  due  to  a  point  f  tree  in  a  semi-infinite  plate  is  used  to  evaluate  the  stresses 
at  the  mid-point  of  each  interval.  Tlie  boi  dry  conditions  that  the  mid-point  of  each  interval  be  free 
from  stresses  are  applied  and  the  values  of  the  body  force  densities  can  then  be  evaluated. 

The  stresses  at  an  arbitrary  point  in  the  plate  and  the  stress  concentration  factor  (SCF)  at 
the  top  or  bottom  of  the  semi-ellipse  or  ellipse  can  be  calculated  by  taking  a  linear  combination  of 
the  body  force  densities. 

Due  to  the  symmetry  of  the  ellipse  or  semi-ellipse  about  the  y-axis,  only  half  of  the 
influence  coefficients  need  to  be  calculated  although,  as  before,  they  must  be  evaluated  for  body 
forces  around  the  total  ellipse. 

3.4.2  Solution  Sequence 

(1)  The  ellipse  or  semi-ellipse  in  a  semi-infinite  plate  is  divided  into  MM  intervals.  For  the 
ellipse,  the  values  of  the  elliptic  parameter  ^  at  the  end  points  of  the  N-th  interval  are 


*N_1  =  -  ff/2  +  2(N-1)tt/MM 


<pN  =  -  tt/2  +  2N  ff/MM 


The  x  and  y  values  of  the  mid-point  of  the  N-th  interval  are 

x  =  acos(-7f/2  +  (2N-l)rr/MM) 
y  =  b  sin  (-ff/2  +  (2N- l)7r/MM) 

For  the  semi-ellipse,  the  values  of  at  the  end  points  of  the  N-th  interval  are 

^N-i  =  (N-  l)rr/MM 
i pN  =  Nct/MM 

The  x  and  y  values  of  the  mid-point  of  the  N-th  interval  are 

x  =  a  cos  ((N-0.5)ir/MM) 
y  =  bsin  ((N-0.5)tt/MM) 

(2)  The  influence  coefficients  </*$.  ,  ayM>  ayM»  rxyM»  rxyNM,  which  are  the  stresses  at  the 

mid-point  of  the  M-th  interval  due  to  a  body  force  on  the  N-th  interval  Y  in  a  semi-infinite 
plate  having  density  px  =  1  or  py  =  1  are  given  in  Appendix  3. 

(3)  The  remainder  of  the  solution  sequence  is  identical  to  that  for  an  ellipse  in  an  infinite  plate. 

3.4.3  Influence  of  the  Integration  Parameter 

The  method  of  numerical  integration  is  identical  to  that  for  an  ellipse  in  an  infinite  plate. 
However,  for  a  semi-infinite  plate  the  influence  coefficients  which  involve  singularities  and  hence 
the  stress  concentration  factors  are  affected  by  the  choice  of  the  value  of  e  in  the  numerical 
integrations. 

Table  VI  shows  the  stress  concentration  factors  obtained  for  various  values  of  the 
integration  parameter  e  for  a  semi-ellipse  with  b/a  =  2.0  in  a  semi-infinite  plate. 

TABLE  VI 


EFFECT  OF  e  ON  SCF  FOR  SEMI-ELLIPSE 
IN  SEMI-IN  FINITE  PLATE 


FRAC 

SCF 

0.010 

5.187062 

0.008 

5.185324 

0.005 

5.182725 

0.001 

5.179275 

extrapolated  to  0 

5.1784 

T 


In  Figure  11  the  SCF  is  plotted  against  e.  The  error  in  the  SCF  is  almost  linear  in  e,  and 
tends  to  5.1784  as  e  tends  to  zero. 

Because  of  the  linearity  of  the  error  in  the  SCF,  it  is  only  necessary  to  obtain  two  or  three 
values  of  the  SCF  for  different  values  of  e  and  extrapolate  to  obtain  the  value  of  SCF  for  e  =  0. 
Another  method  of  approximating  the  true  SCF  for  a  given  number  of  cuts  is  to  choose  e0  small 
enough  that  the  difference  between  the  value  of  the  SCF  at  e  =  0  and  the  value  at  e0  is  negligible. 

3.4.4  Results 

Tables  VII  and  VIII  compare  the  results  obtained  for  a  semicircular  notch,  and  an  ellipse 
with  b/a  =  2.0  and  e  =  0.5,  with  those  of  Nisitanil 1  for  several  values  of  the  number  of  intervals 
MM.  The  agreement  is  excellent.  Discrepancies  can  be  explained  by  the  fact  that  the  integration 
parameter  e  is  taken  to  be  very  small,  but  no  extrapolation  has  been  done  to  determine  the  true 
value  at  e  =  0.  This  was  done  to  reduce  the  amount  of  computer  time  used. 

TABLE  VII 


STRESS  CONCENTRATION  FACTORS  FOR  A 
SEMICIRCLE  IN  A  SEMI-INFINITE  PLATE 


'^TB 

MM 

Present  Results 

Nisitanil4' 

16 

3.053872 

3.0527 

24 

3.057044 

3.05595 

32 

3.059274 

3.05908 

48 

3.061416 

3.06123 

a  =  1;  b  =  1;  FRAC  =  0.001 


TABLE  VIII 


STRESS  CONCENTRATION  FACTORS  FOR  AN 
ELLIPSE  IN  A  SEMI-INFINITE  PLATE 


MM 

Present  Results 

Nisitanil  *  1 

16 

5.696703 

— 

24 

5.715367 

5.7147 

32 

5.724825 

5.7243 

48 

5.734337 

5.7340 

a  -  1;  b  =  2;  e  =  0.5;  FRAC  =  0.001 


Figures  12  and  13  show  KTA  and  KTB  as  a  function  of  1/MM.  The  error  in  the  KT’s  is 
almost  linear  in  1/MM.  For  the  semicircle  on  the  edge  of  the  plate,  KTB  tends  to  3.0650  as  MM  tends 
to  infinity,  compared  with  a  value  of  3.0654  obtained  by  Nisitani*4'  and  3.0653  by  Isida^5 1 .  For  the 
ellipse  (b/a  =  2.0,  e  =  0.5)  KTA  tends  to  5.7531  as  MM  tends  to  °°,  compared  with  a  value  of  5.753 
obtained  by  Nisitani. 


3.5  Combination  of  a  Crack  and  an  Elliptical  Hole  in  a  Semi-infinite  Plate 
3.5.1  Formulation 

Figure  14  shows  the  case  of  a  crack  and  an  elliptic  edge  notch  in  a  semi-infinite  plate.  As 
before,  the  semi-ellipse  is  divided  into  Ml  equal  intervals  and  the  crack  into  M2  equal  intervals,  with 
the  body  force  densities  set  constant  within  each  interval.  The  stress  field  due  to  a  point  force  in  a 
semi-infinite  plate  and  its  derivatives  are  used  to  calculate  the  stresses  at  the  mid-point  of  each  interval 
on  both  the  crack  and  the  semi-ellipse  due  to  the  body  forces  on  both,  and  the  conditions  to  make 
the  mid-point  of  each  interval  free  of  stress  are  applied.  The  body  force  densities  can  then  be 
determined  as  before. 


3.5.2  Solution  Sequence 

(1)  The  semi-ellipse  is  divided  into  Ml  equal  intervals.  The  end  points  of  the  N-th  interval 
have  i p  values 


*N-1  =  (N-1)tt/M1 


«  Nir/Ml 


where  i p  is  the  elliptic  parameter  such  that 

£  =  a  cos  i p;  r/  =  b  sin  y 

The  x  and  y  values  of  the  mid-point  of  the  M-th  interval  are 

x  =  a  cos  ((M- 0.5) 7r /Ml) 
y  =  b  sin  ((M- 0.5)  jt/MI) 

The  crack  is  divided  into  M2  equal  intervals.  The  end  points  of  the  N-th  interval  have  values 
on  the  y  axis 

r?N_!  =  e-c  +  2(N-l)c/M2 


rjN  =  e-c  +  2Nc/M2 


where  e  is  the  distance  from  the  edge  of  the  plate  to  the  centre  of  the  crack,  and  c  is  half 
the  length  of  the  crack. 
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The  mid-point  of  the  M-th  interval  has  a  value 

y  =  e-c  +  (2M-l)c/M2 


(2)  The  influence  coefficients 


XN  XN  YN  YN  ^XN 
°xMe ’  °yMe ’  uxMe>  uyMe 


xyMe  * 


r 


YN 

xyMe 


a 


yN 

xMe’ 


ryN  yN 
yMe*  'xyMe 


XN  YN 
xMc»  °xMc 


and 


a 


yN 

xMc 


are  given  in  the  appendices. 

The  first  six  coefficients  are  the  stresses  at  the  mid-point  of  the  M-th  interval  on  the  ellipse 
due  to  point  forces  on  the  N-th  interval  of  the  ellipse  having  densities  px  -  1  or  py  =  1 
(see  Appendix  3). 

The  seventh  to  ninth  coefficients  are  the  stresses  at  the  mid-point  of  the  M-th  interval  on 
the  ellipse  due  to  a  pair  of  point  forces  in  the  N-th  interval  of  the  crack  having  density 
7  =  1  (see  Appendix  6). 

Hie  10th  and  11th  coefficients  are  the  stresses  at  the  mid-point  of  the  M-th  interval  on  the 
crack  due  to  point  forces  on  the  N-th  interval  of  the  ellipse  (see  Appendices  3  and  6). 

Hie  12th  coefficient  is  the  stress  at  the  mid-point  of  the  M-th  interval  on  the  crack  due  to 
a  pair  of  point  forces  on  the  N-th  interval  of  the  crack  (see  Appendix  2). 

(3)  The  boundary  conditions  that  the  mid -point  of  each  interval  on  the  ellipse  and  crack  be 
free  from  stress  are  applied: 


Ml 

E 

N=  1 


M2 

Tn  (°xMe  COS  6  +  TxyNMe  sin  6)  +  <  COS  0  =  0 

N  =  1 
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M1 

^  [pxN  (°yMe  sin  0  +  TxyMe  COS  6)  +  ^  (°yMe  sin  0  +  'xyMe  COS  #)] 


23  7~  (‘ 


<£2.  ™  9  ♦  r*yMe  COS  * 


’)  '0 


Ml 

£[- 


„  XN  .  „  _YN 
PxN  °xMc  +  PyN  ffxMc 


M2 

■]*E 


°lL  +  °l 


The  result  is  a  system  of  (2M1 
and  y. 

3.5.3  Numerical  Integration 


+  M2)  linear  equations  in  (2M1  +  M2)  unknowns  px, 


In  the  numerical  integration  of  the  influence  coefficients,  singularities  occurred  in  the 
evaluation  of  the  first  six  coefficients  discussed  in  3.5.2,  and  in  the  12th.  The  procedure  for  handling 
the  first  six  coefficients  is  discussed  in  3.3.3,  while  the  procedure  for  the  12th  coefficient  is  discussed 
in  3.1.5.  However,  it  is  not  necessary  to  evaluate  the  12th  coefficient  numerically.  As  seen  in  3.1.3 
(for  a  crack  of  length  2c)  or  3.2.1  (for  a  crack  at  the  edge  of  a  semi-ellipse)  the  12th  coefficient  can 
be  evaluated  analytically. 

3.5.4  Results 

The  values  of  the  correction  factors  at  both  tips  of  a  crack  in  a  semi-infinite  plate  having  a 
semi-elliptical  notch  are  given  for  various  values  of  the  number  of  intervals  Ml  and  M2  in  Table  IX. 


TABLE  IX 

CORRECTION  FACTORS  FOR  A  CRACK  NEXT  TO  A 
SEMI-ELLIPTICAL  NOTCH  IN  SEMI-INFINITE  PLATE 


No.  of  Intervals 


Correction  Factors 


Ml 

M2 

Fia 

F.b 

16 

16 

2.201207 

1.935246 

24 

24 

2.211444 

1.938788 

32 

32 

2.217546 

1.940605 

48 

48 

2.222124 

1.941886 

c  =  100.0;  p  =  1000.0;  c/p  =  0.1 
a  =  1000.0;  b  =  1000.0;  b/a  =  1.0 
c/h  =  0.4;  FRAC  =  0.05;  e  =  FRAC  •  tt/MM 


w; 


v  I-VvV 
-  •  ^  =  - 


.  -  .v 


The  correction  factors  are  plotted  against  1/(M1  +  M2)  in  Figure  15  and  again  the  error  is  almost 
linear  in  1/(M1  +  M2).  The  extrapolated  values  as  Ml  +  M2  tend  to  infinity  are  2.2329  for  FIA  and 
1.9454  for  F1B  compared  with  the  values  of  2.232  and  1.944  obtained  by  Nisitani  l 1 1 . 

4.0  SUMMARY  AND  CONCLUSIONS 

An  outline  of  Nisitani’s  body  force  method  has  been  presented  and  several  problems  have 
been  solved  using  the  technique.  The  problems  include  single  embedded  cracks  and  edge  cracks  and 
single  elliptical  and  semi -elliptical  notches  in  semi-infinite  plates,  as  well  as  the  problem  of  an 
embedded  crack  near  a  semi-elliptical  notch  in  a  semi-infinite  plate.  Results  agree  very  well  with  those 
of  Nisitani  and  solution  times  are  relatively  short. 

Details  of  the  numerical  technique  have  been  presented,  including  a  brief  discussion  of  the 
computer  precision  required  and  the  method  for  integrating  over  the  singularities.  In  addition,  all  the 
analytical  details  for  each  problem  have  been  presented,  gathering  together  information  from  several 
of  Nisitani’s  papers,  many  of  which  are  in  Japanese. 

Full  details  of  the  computer  programs  are  provided  in  a  separate  report. 
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FIG.  7:  CRACK  IN  SEMI-INFINITE  PLATE 
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FIG.  9:  EDGE  CRACK  IN  SEMI-INFINITE  PLATE 
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FIG.  12:  ELLIPSE  IN  SEMI-INFINITE  PLATE 
(b/a  *  2.0,  e  =  4.0,  FR AC  *»  0.001  > 
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FIG.  13:  SEMI-ELLIPSE  IN  SEMI-INFINITE  PLATE 
(b/a  -  1.0,  a  -  0,  FRAC  -  0.001) 
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APPENDIX  1:  FUNDAMENTAL  STRESS  HELDS 


For  plane  stress,  the  stress  fields  at  (x,y)  for  point  forces  X  and  Y  acting  at  (£,tj)  are  l 


o*  =  -  FC  {(3  +  v)i2  +  (1  -  y)m2}  X 
o*  =  F £  {(1  -  v)£2  -  (1  +  3t>)m2jx 
rfy  =  -  Fm  {(3  +  v)i2  +  (1  -  i>)m2}x 


o*  =  -  Fmj(l  +  3k)£2  -  (1  -^m^Y 
aj  *  -  Fm  1(1  -  v)i2  +  (3  +  i>)m2}  Y 
Txy  “  “  FC  ((*  "  v)*2  +  (3  +  t>)m2  J  Y 


for  an  infinite  plate,  and 


of  =  [-  FC  ((3  +  v)9.2  +  (1  -  i>)m2|  -  GjC  ((3  +  v)&2  +  (1  -  »»)n2} 

-  G2  {(1  -  v)  (C5  +  2  (n2  -  n)C3  +  n2  (n2  -  2n)C) 

+  (1  +  v)  (-  2  (n2  -  1)C3  +  n2  (-  2n2  +  8n  -  6)£)}]  X 

Oy  =  [f£  {(1  -  v)i2  -  (1  +  3y)m2 1  +  GjC  ((1  -  v)i2  -  (1  +  3v)n2} 

-  G2  ((1  -  v)  (C5  +  2  (n2  +  n)C3  +  n2  (n2  +  2n)C) 

+  (1  +  v)  (2  (-  n2  +  2n  -  1)£3  +  n2  (-  2n2  -  4n  +  6)G)(J  X 

rfy  *  [-  Fm  ((3  +  v)i2  +  (1  -  i>)m2}  -  G,n  {(3  +  v)l2  +  (1  -  i>)n2} 

+  G2  {(1  -  v)(-  E4  +n4)  +  (1  +i>)(6n(n-l)£2  +  n3  (-2n  +  2))(]  X 


o*  =  [-  Fm  {(1  +  3p)C2  -  (1  -  p)m2|  +  Gjn  |(1  +  3p)£2  -  (1  -  v)n2} 
-G2((l  -  i>)((n  +  l)?4  +  2n3 £2  +  n3  (n2  -  n)) 


+  (1  +  v)  (2C4  +  n  (2n2  +  6a  -  6)£2  +  n3  (2n2  -  4n  +  2))|]  Y 


a*  =  Fm  j(l  -  v)£2  +  (3  +  i»)m2}  +  G^  {(1  -  v)£2  +  (3  +  i>)n2  J 
-  G2  {(1  -  v)  ((n  -  lJC4  +  2n3  £2  +  n3  (n2  +  n)) 

+  (1  +  v)  (n  (2n2  -  6n  +  6)£2  +  n3  (2n2  +  2n  -  2))|]  Y 


Txy  =  ["  FC  {(1  '  V)S2  *  <3  +  +  Gi fi  ((1 '  +  (3  +  "fr2) 

-  G2  ((1  -  v)  (2n£3  +  2n3  £)  +  (!+  v)  (2£3  +  n2  (8n  -  6)£))]  Y 


for  a  semi-infinite  plate,  where 


£  = 


m 


y  -  v 


n  = 


y +  n 


4jry(fi2+m2)2  1  47ry(£2+n2)2 


G2  " 


2wy  (£2  +n2)3 


APPENDIX  2:  ANALYTIC  INTEGRATION  FOR  CRACK  IN  SEMI-INFINITE  PLATE 


Stress  at  the  M-th  interval  due  to  a  pair  of  point  forces  at  the  N-th  interval  on  a  crack: 


,iN  = 

JxM 


N  ^  |X=1 

N  I x=f=0 


{4>n  Vb2  -  (i?  -  e)2  ) 


dr) 


>N  =  1 


f N  I  1  J_  I  1  12y  12y2 

J  4ir  (y  -  n)2  4M  (y  +  n)2  (y  +  n)3  (y  +  t»)4 

^N-l  *“ 


X  {4  Vb2  -  ( r)  -e)2|  dr) 


where  tjn-1  and  tjn  denote  the  values  of  tj  at  both  ends  of  the  N-th  interval.  Die  integral  can  be 
obtained  in  closed  form  as  follows^ 1 ' : 


°xM  =  p<%>  -f(t?n-i) 


in  which  F(t?)  stands  for 

1  K/b2  -  (tj  -  e)2  .  17-e  y-e 

*|_  y-7*  b  Vb2  -  (y  -  e)2 


X  In 


b2  -  (y~e)(7?  -e)  +  V{b2  -  (y-e)2}{b2  -  (r?-e)2| 
b(y-77) 


l[  Vb2  -  (r)  -  e)2  .  t?  -  e  y  +  e  .  b2  +  (y +  e)(r?  -  e)l 

w  L  y  +  1?  b  V(y  +  e)2  -  b2  b  (y  +  tj)  J 

+  6y  f  Vb2  -  (tj-e)2  f  _  (y  +  e)(y +  tj)|  +  b2  x  gin_i  b2  +  (y  +  e)  (r)  -  e)l 

L  (y  +  t?)2  I  (y  +  e)2-b2  I  {(y  +  e)2  -  b2}*  *”  b(y  +  T,)  J 

_  2y^  |>/b2  -  (1?  -  e)2  x  I  2  +  (y +  e)(y +  r;)  +  {2b2  +  (y  +  e)2}  (y  +  tj)2) 
w  [_  (y  +  T?)3  I  (y  +  e)2  -  b2  I(y  +  e)2-b2}2  f 

+  3b2  (y  +  e) 

{(y  +  e)2  -b2}» 


sin 


-1  b2  +  (y +  e)(r?- 

b(y  +  17) 


e) 


w//.' 
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APPENDIX  3:  INFLUENCE  COEFFICIENTS  FOR  ELLIPSE 


I 


i 


R 


£ 


i 


The  influence  coefficients  can  be  calculated  from  the  stress  fields  for  a  point  force  X  or  Y 
acting  at  (£,jj)  (Appendix  1).  Since  the  relations  drj  =  b  cos  d and  -  d£  =  a  sin  dyJ  prevail  along 
the  ellipse,  the  coefficients  can  be  expressed  as^1 ' : 

°xM  =  /  (S.i7,x,y)jx=1  bc°s^d^ 

N 

°yM  “  J  °y  <^,T?,X,y)jx=  1  bCOS^d^ 

N 


M  =  /  (^.X.y)|x=i  bcos^d^ 

N 

and 

°xm  =  f  °*  (*.n,*,y)|v=1asin*dv> 
N 


ayM  =  J  ay  (^.»?.X»y)|Y=l  »sin^d^ 

N 

TxyM  =  /  Txy  (U,*,y)|Y=  1  8  Sin  * 

N 

where  /N  stands  for  integration  of  the  N-th  interval  and  a*,  a*,  o*.  ay,  rfy ,  r*y  are  given  in 
Appendix  1. 


i 


APPENDIX  4:  INTEGRALS  OVER  SINGULARITIES  FOR  ELLIPSE 


The  integrals  over  the  singularities  for  the  ellipse  can  be  approximated  by  the  following 
simple  analytic  expressions^6) : 

Aa™  =  e - - - [6w6  +  (3  -  2k2)  w4  -  6k2w2  -  k4] 

4ir(k2+W2)3 


AaXM  =  e - [-  2w6  -  (1  -  6k2) w4  +  6k2w2  -  k4] 

yM  4ir(k2  +  w2)3 


Ar™  =  e - [-9W5  -  (6  +  2k2)w3 +k2  (2  -  k2)w] 

4jr(k2  +  w2)3 


Aa™  =  e - [~  w6  +  6k2 w4  +  k2  (6  -  k2)w2  -  2k4] 

4ff(k2  +  w2)3 


where 


Aa™  =  € - [-  w6  -  6k2w4  +  k2  (2  -  k2)w2  +  6k4] 

4ir(k2  +  w2)3 

AtxvM  =  ~€ - - - [(1  -  2k2)  w5  +  2k2  (1  +  3k2)  w3  +  9k4w] 

4jr(k2  +  w2)3 

k  =  b/a 

w  =  tan  <fi0  (<fi0  =  elliptic  parameter  value  of  mid-point  of  M-th  interval) 


Note  that  there  is  a  sign  change  from  the  original  (Reference  6)  in  the  expression  for  A  a 


APPENDIX  5  :  STRESSES  CLOSE  TO  ELLIPSE 


The  stresses  at  a  point  infinitesimally  close  to  the  ellipse  subjected  to  a  body  force  may  be 
calculated  in  the  following  manner  as  given  by  Nisitanil4  *  in  somewhat  less  detail. 

Consider  an  elliptic  arc  (Fig.  Al)  of  length  2e  on  the  x'  axis.  When  a  body  force  having 
densities  py«,  px-  acts  on  the  arc,  the  stresses  at  a  point  (0,b)  (b  small)  due  to  px  dx'  and  py  dx'  are 
obtained  by  differentiating  the  equations  of  the  stress  field  in  an  infinite  plate  (Appendix  1). 


d<V 


x'  (3x'2  +  b2) 
4ff(b2  +  x'2)2 


+ 


b  (b2  -  x'2) 
4tt  (b2  +  x'2) 


Pydx' 


da 


y' 


x'(x’2-b2)  ,  b  (3b2  +  x'2) 

- px,dx - 

4jr(b2  +  x'2)2  4ff(b2+x'2)2 


dr 


x  y 


b(3x'2  +  b2) 
(b2  +  x'2)2 


px-dx' 


x  (3b2  +  x'2) 
+  — - - - 

(b2  +  x'2)2 


Py-dx' 


Putting  x'  =  b  tan  0 ;  =  tan-1  — ;  dx'  =  b  sec2  0  Ad  and  integrating: 


ip  r- 

1  C  .  .  n  (z  sin2  6  \  ,  , 

=  —  I  sm  6  cos  6  I - +  l|  px-  +  (cos2  9  -  sin2  6)  py>  Ad 

4  L  \  c°s2  6  ) 


i  v 

ay>  =  —  J*  |^-  sin  9  cos  9  (tan2  0-1)  px.  -  (3  cos2  9  +  sin2  0)  py>J 


d0 


1  ^ 

•y‘  -  —  j  (3  sin2  0  +  cos2  0)  px«  +  sin  0  cos  0  (3  +  tan2  0)  py-J 


d0 


-<P 


As  b  tends  to  zero,  0  tends  to  tt/2  and  the  stresses  become: 


ox-  *  0;  av' 


1 

2 


"  77  Py'  I  Tx'y' 


1 

2 


(Al 


Note  the  difference  in  sign  of  r.>v>  compared  with  that  in  Reference  4. 


>  V  V 


Resolving  the  force  densities  along  x  and  y  axes  as  in  Figure  A2  gives 


px-  =  -  py  sin  0  cos  0  +  px  cos  0  sin  0 


-(px  -  py )  sin  20 


py'  =  Py  sin2  0  +  px  cos2  0 


Resolving  the  stresses  along  the  x'  and  y'  axes  as  in  Figure  A3  gives 


A ox  =  ox-  sin2  0  +  oy'  cos2  0  +  2rxy  sin  0  cos  0 


Aoy  *  oX'  cos2  0  +  oy'  sin2  0  -  2rxy  sin  0  cos  0 


Ar  =  -  ox.  cos  0  sin  0  +  oy-  cos  0  sin  0  -  Txy  (cos2  0  -  sin2  0) 


Substituting  Equations  (Al)  and  (A2)  into  the  above  gives 


_  1  ,  .  ,  „  ,  sin  0  cos  0 

Aox  *  -  —  cos-4  0  (px  cos*  0  +  py  sin*  0) - - - (px  -  p  )  sin  20 


-  —  px  cos2  0  (1  +  sin2  0)  +  ^  py  (cos2  0  sin2  0) 

«  A 


Px  Py 

—  (5  +  4  cos  20  -  cos  40)  +  —  (1  -  cos  40) 


a  1  .  j  ,  ,  .  .  ,  sin  0  cos  0 

Aoy  -  -  —  sin*  0  (px  cos*  0  +  py  sin*  0)  + - - - (px  -  py)  sin  20 


px  (cos2  0  sin2  ~  Py  (sin2  0)  (1  +  cos2  0) 

a  2 


Px  Py 

—  (1  -  cos  40)  -  —  (5  -  4  cos  20  -  cos  40 ) 


APPENDIX  6:  ADDITIONAL  INFLUENCE  COEFFICIENTS  FOR 
COMBINED  CRACK  AND  ELLIPSE 


The  stresses,  o^e,  o^e  and  7^”Me,  on  t^le  M-th  *nterval  °f  30  ellipse  due  to  a  pair  of 
point  forces  on  the  N-th  interval  of  a  crack  are  given  as^7' : 


-s 


**  t-o 


'c2  —  (17  —  e)2  dr? 


,7N  = 


J  at 


a*  c-o 


'c2  —  (17  —  e)2  dr? 


TrN  _  I  — 
*xyMe  J  3t 


'c2  -  (r?  -  e)2  dr? 


where 


— =  - [_  3A4  +  6A2ro2  +  m4]  + - 

t=o  4rry2  (A2  +  m2)3  4rry2  ( A2  +  n2 )4 


X  [-  5A6  +  (13n2  +  12n  -  12)  A4  +  n4  (17n4  -  72n  +  72)  A4  -  n4  (n4  -  12n  +  12)] 


^4  4.  _2  i-in~2 


2  -  r,4  r«2  _ 


9ffv  1  1 

-f  =  - [A4  -  6A2m2  +  m4]  + - 

t=o  4rry2  (A2  +  m2)3  4rry2  (A2  +  n2)4 

X  [-  A6  +  (5n2  -  36n  +  12)  A4  +  n2  (5n2  +  56n  -  72)  A2  -  n4  (n2  +  4n  -  12)] 


®  xy  Am  9  ~  1 

— —  =  - [-  3A2  +  m2]  + - 

^  i=o  2rry2  (A2  +  m2)3  4rry2  (A2  +  n2)4 

X  [2  (3n  +  2)  A5  +  4n  (n2  -  14n  +  12)  A3  -  2n3  (n2  -  18n  +  24)  A] 


x  v-y  v+y 

A  =  — ,  m  =  - ,  n  *  - 

y  y  y 


ar* 

xy 

Note  the  change  of  sign  in  — —  as  compared  to  the  equation  given  in  Reference  7. 
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